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ABSTRACT 

Tomography of redshiftcd 21 cm transition from neutral Hydrogen using Fourier synthesis telescopes 
is a promising tool to study the Epoch of Reionization (EoR). Limiting the confusion from Galactic 
and Extragalactic foregrounds is critical to the success of these telescopes. Instrumental response or 
the Point Spread Function (PSF) of such telescopes is inherently 3 dimensional with frequency mapping 
to the Line of Sight (LOS) distance. EoR signals will necessarily have to be detected in data where 
continuum confusion persists; therefore, it is important that the PSF has acceptable frequency structure 
so that the residual foreground does not confuse the EoR signature. This paper aims to understand the 
3 dimensional PSF and foreground contamination in the same framework. We develop a formalism to 
estimate the foreground contamination along frequency, or equivalently LOS dimension, and establish a 
relationship between foreground contamination in the image plane and visibility weights on the Fourier 
plane. We identify two dominant sources of LOS foreground contamination — 'PSF contamination' and 
'gridding contamination'. We show that 'PSF contamination' is localized in LOS wavenumber space, 
beyond which there potentially exists an 'EoR window' with negligible foreground contamination where 
we may focus our efforts to detect EoR. 'PSF contamination' in this window may be substantially reduced 
by judicious choice of a frequency window function. Gridding and imaging algorithms create additional 
'gridding contamination' and we propose a new imaging algorithm using the Chirp Z Transform (CZT) 
that significantly reduces this contamination. Finally, we demonstrate the analytical relationships and 
the merit of the new imaging algorithm for the case of imaging with the Murchison Widefield Array 
(MWA). 

Subject headings: Techniques: interferometers — methods: data analysis — methods: analytical — cosmology: ob- 
servations 



1. Introduction 

Several synthesis telescopes such as MWA0, LO- 
FAR@, and PAPEE| are currently being built to study 
the Epoch of Reionization (EoR) by observing the red- 
shifted 21 cm line transition from neutral hydrogen. 
While future instruments may map the 21 cm bright- 
ness in three dimensions (frequency corresponds to 
line of sight), current experiments lack the required 
sensitivity and instead hope to estimate the power 
spectrum of 21cm brightness fluctuations. The bright- 
ness of these fluctuations is about 5 orders of mag- 
nitude smaller than the Galactic and Extragalactic 
foregrounds. Moreover, removing the effects of the 
telescope point spread function (PSF) from the syn- 
thesis images may not be possible to the level of sen- 
sitivity demanded by EoR power spectrum measure- 
ments. Consequently, foreground contamination and 
instrumental effects present a major hurdle to 21 cm 
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tomography. 

A host of simulations have demonstrated the poten- 
tial in current instruments to recover the EoR power 
spectrum despite thermal uncertainties, foreground 
contaminants, and instrumental effects. Among in- 
strumental effects, simulations have shown the fre- 
quency dependence of PSF to be o f particular concern 
in EoR power spectrum estimation (jBowman. Morales fc Hewitt 
2009). In 21 cm tomography, frequency maps to the 
line of sight, and any instrumental variation in fre- 
quency can potentially 'mimic' line of sight 21 cm 
brightness fluctuations. Most of the proposed fore- 
ground subtraction algorithms rely on the ability 
to distinguish smooth spectral features of Galactic 
and Extragalactic contaminants from the relatively 
rapid spectral fluctuations (expected from theory) 
of the 21 cm EoR signal. Spectral features in the 
telescope PSF may affect our ability to 'fit' and 
remove undesirable foregrounds. The sidelobe re- 
sponse of confusing foreground contaminants may 
vary with frequency and generate an instrumen- 
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tal spectral structure at each image pixel — a phe- 
nome non that has been dubbed 'mode-m ixing' in liter- 
ature ( Bowman. Morales fc Hewitt 2009). Apart from 
mode-mixing, calibration errors and choice of gridding 
algorithms can also lead to undesirable instrumental 
response in frequency domain. 

While simulations have evaluated such mode- 



mixing effects for realistic arrays like MW A ([Liu et al, 
20091: iBowman. Morales fc Hewitt! l2009t) . this paper 
presents an approach that is analytical, and hence, 
more general. In particular, the relationship between 
the PSF in image coordinates (Imf) and the telescope 
sampling function in Fourier coordinates (uvf) is es- 
tablished, and the mode-mixing problem is cast in 
Fourier space where the telescope measurements are 
really made. In this formalism, implications for choice 
of experimental parameters and foreground subtrac- 
tion algorithms are discussed. Apart from the an- 
alytical treatment, various practical data processing 
operations that influence the PSF are discussed using 
an example of snapshot imaging with MWA. 

The remainder of the paper is organized as follows. 
Section [5] defines the problem analytically and derives 
an expression relating the line of sight foreground con- 
tamination to the weighting function in Fourier space. 
To better understand the relationships established in 
Section [21 Section [3] applies them to the simple case of 
imaging a point source with a complete visibility dis- 
tribution. The analytical treatment in Section [5] does 
not take into account the effects of discretization of 
visibilities in the re-gridding and imaging algorithms. 
Section |4] discusses these effects and complements the 
results of Section [21 In particular, Section [4] describes 
a new imaging algorithm based on the Chirp Z Trans- 
form (CZT) that alleviates the problem of foreground 
contamination. Section [5] demonstrates the concepts 
developed in Sections [21 02 and [4] through simulating a 
snapshot observation with MWA. Section [6] draws up 
conclusions and proposals for future work. 

2. Analytical treatment 

We first analytically evaluate the effect of a fre- 
quency variant PSF at a given image pixel due to an 
unsubtracted point source by evaluating its contribu- 
tion to foreground contamination in LOS wavenum- 
ber space. We call this contamination 'PSF contam- 
ination'. A single point source is considered here to 
attain an intuitive understanding of 'PSF contamina- 
tion'. Generalization for a realistic sky will be made 
in subsequent sections. Array element locations are 
assumed to be coplanar. 

2.1. Notation 

We define the co-ordinate system to avoid ambigu- 
ity. The baselines are defined on a local plane whose 



origin is at the latitude and longitude of the array cen- 
ter. A two dimensional vector x = [x y] (in meter 
units) represents the projected baseline on the local 
plane. The corresponding baseline u = [u v] (in wave- 
lengths) at a frequency / is then given by u = x.(//c), 
where c is the speed of light. A two-dimensional direc- 
tion vector 1 towards any direction is represented by 
its direction cosines: 1 = [I to]. An interferometer ob- 
servation gives a visibility cube with two spatial axes 
representing the baseline and a frequency axis. An 
image cube on the other hand has two axes represent- 
ing direction cosines on the sky and a third frequency 
axis. Assuming co-planarity of array elements, within 
the framework of Fourier synthesis imaging, u and 1 
are Fourier conjugates. 

2.2. Weighting function 

The telescope sampling function is represented as 
W(x, /), which is the weight assigned to the baseline 
represented by the vector x, and includes any visibility 
taper applied prior to imaging. In general, W(x) de- 
pends on instrumental, observational and processing 
parameters. Among instrumental parameters is the 
antenna layout that determines the set of available 
baselines. Observational parameters, which include 
pointing direction, duration of rotation synthesis, and 
mode of observation — drift scan or track mode, define 
the uv coverage. Observational parameters along with 
the antenna layout determines the 'natural' weight at 
each x — number of visibilities falling into the pixel cen- 
tered at x. Among processing parameters are the vis- 
ibility grid size and extent, convolution kernel used to 
grid the visibilities, and any visibility mask and/or ta- 
per used during imaging. The weighting function may 
be factored into a frequency independent geometric 
component, a frequency dependent 'taper' component, 
and a frequency window function: 



VF(x,/) = VF(x)T(x,/)i?(/). 



(1) 



Instrumental and observational parameters are geo- 
metric in nature and hence frequency independent. 
They may be sufficiently represented by W(x). Pro- 
cessing parameters like visibility 'taper' function and 
visibility masks are frequency dependent and are rep- 
resented by T(x, /). B(f) is the window function in 
frequency assigned to each visibility. 

2.3. Derivation 

An image made with a weighted baseline distribu- 
tion W(x 7 f) has a frequency dependent PSF A (1, /), 
where 1 is the two dimensional unit vector representing 
direction on the sky. Tilde notation is used for image 
domain functions as opposed to their Fourier domain 
counterparts. Consider the sidelobe contribution at 
zenith (1 = [0 0]) from a source with flux density 
of unity and spectral index of zero located at any 1. 
This contribution 5(1, /) varies with frequency due to 
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a changing PSF. We wish to compute the energy in 
these spectral variations on different scales in terms of 
the telescope weighting function W(x., /). 

PSF A(l, /) is the Fourier Transform of the uv- 
weighting function A(u, /). We may represent this 
as 

f°° du f°° dv A(u, f) cxp(— i27ru.l) 

{U) j\duf\d,,i,,.n ' {Z) 

where the denominator is included to obtain A(l, f) as 
the normalized PSF for all frequencies, and essentially 
does bandpass calibration of visibilities. The inclu- 
sion of A(l, /) essentially does bandpass calibration of 
visibilities and results in a PSF that is normalised at 
all frequencies. Any additional frequency weighting 
may be absorbed into the frequency window function 
B(f). The function A is just a scaled version of the 
function W through the scaling relationship given by 
u = x.(//c). We may represent this scaling relation- 
ship as 



A(u,f) = W\uj,f) =W(x,/). (3) 

The quantity uc/f has units of distance and is fre- 
quency invariant: it is simply the baseline distance 
in meters. We use this to make the substitution 
x = uc/f, y — vc/f and get 
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where Aw = dx j_ OQ dy W(x, /) is the area un- 
der W. The importance of this substitution is that it 
expresses all baselines in units of meter as opposed to 
wavelength and makes them independent of frequency. 

Consider a point source with a flux density of unity 
and a spectral index of zero located at 1. The response 
at zenith due to a sidelobe on this source is given by 



5(1, /) = !(!,/) 



(5) 



We are interested in the energy in 5(1, /) at different 
frequency scales. In other words, to evaluate 'PSF con- 
tamination', we have to compute the power spectrum 
5(1,77) by Fourier transforming 5(1,/) with frequency 
offset A/ and time 77 as Fourier conjugates. We use 
A/ instead of / because the spatial EoR fluctuations 
are evaluated about an origin tq corresponding to a 
frequency /o and redshift zq. 

5(77, 1) may be expressed as 

y* 00 

5(1,77)= / d(A/)B(/)i(l,/)ex P H2^7A/). (6) 



Using A/ = / - f a , f = (f L + f H )/2, where f L and 
fn are the lowest and highest frequencies in the instru- 
ment bandpass, and writing the conjugate variable as 
/ instead of A/ we get 

/oo 
df 
-OO 

£(/)i(l,/)expH2^7/). (7) 

Substituting from Equation U for A(l, f) and inter- 
changing the order of integrations we get 

/•oo />oo 

5(1,77) =exp(727r77/ ) / dx dyW(x) 



(8) 



#^W)exp fi27r/(— -7, 
A w \ \ c 



The inner integral, in general, may have to be eval- 
uated numerically. However, to have an insight into 
the nature of 5(1, 77), we evaluate Equation [5] for the 
case where there is no taper or mask applied in the 
uv domain. In this case, T(x, /) = 1, and Aw is fre- 
quency independent. The inner integral is now just a 
scaled and shifted version of the Fourier transform of 
the frequency window function. Consider a rectangu- 
lar frequency window defined by 



B(f) 



\ i£f€[f L ,f H ] 
otherwise, 



(9) 



where B = fjj — fL is the frequency bandwidth. The 
inner integral now evaluates to 
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This gives us the relationship we are looking for: 
exp(i27T77/ ) 



5(1,77) 



Aw 



dx 



dy 



(11) 



The power spectrum 5(1, 77) is a convolution of the 
weighting function W(x) and the convolving kernel B, 
which is a shifted version of the Fourier transform of 
the frequency window function B(f). 

3. Ideal case of complete visibility coverage 

This section provides an insight into the terms and 
equations of the analytical derivation, by understand- 
ing them in the case of a visibility distribution with 
complete coverage. Complete coverage refers to the 
case where there exists at least one visibility measure- 
ment in every pixel in the uv grid at all frequencies 
within the instrument bandwidth, and over the spa- 
tial frequency range spanned by the array configura- 
tion. A key idea in the analytical treatment of Section 
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Fig. 2. — Plot depicting the operation in Equations [TOI 
and [11] See Section 13.11 for a detailed explanation of 
the operations involved. 



Fig. 1. — Plots showing an example of PSF sidelobe 
structure in frequency for a uniformly weighted visi- 
bility distribution with complete coverage (W(x) = 1 
for |jc| < 1000, W(x) = otherwise). Top left panel 
shows a small section of the PSF in the vicinity of 
(Z , f a ) = (0.1, 200 MHz). Top right panel shows PSF 
variation at l D as a function of fractional frequency, 
Sf/fo- Bottom left panel shows PSF variation at f„ as 
a function of fractional directional cosine, Sl/l . Bot- 
tom right panel shows the PSF Fourier transformed 
with respect to frequency to reveal the 'PSF contami- 
nation' in ry space due to a point source at l . 

[2] is the invariance of the PSF in 1/. This concept is 
manifested in the invariance of Equation [4] in 1/ for 
a frequency invariant weighting. Frequency invariant 
weighting implies a 'taper' function T(x) that is in- 
variant in frequency and dependent only on x. Figure 
[TJ shows an image of the PSF in If domain for a uni- 
formly weighted visibility distribution with complete 
coverage. The plot extends over a 10% spread in both 
frequency and direction cosine (Al/l Q — Af / f a — 0.1). 
Due to the If invariance, we expect the PSF to have 
the same value at all points defined by If = l f - 
Consequently, the variation in an interval 81 of the 
PSF evaluated at f a , is the same as the variation in 
an interval f a of the PSF evaluated at l , so long as 
Sl/l = —Sf/f or approximately 5l/l = —Sf/f a . This 
can be obtained by differentiating If = l f and then 
approximating / to f a and I to l a , which hold for 
51/1 << 1 and 6f/f << 1. Figure [TJ also evaluates 
'PSF contamination' in the zenith pixel in rj space due 
to PSF sidelobes on a 1 Jy point source at l a . The lo- 
calization of this contamination in light of Equations 
[11)1 and [TT1 is pivotal to our understanding of 'PSF con- 
tamination' and deserves a detailed explanation. 



3.1. Localization of 'PSF contamination' 

Consider a single point source at direction cosine 
l entering the zenith pixel through the PSF sidelobe 
at l a . This is depicted in Figure [2] Due to the 1/ 
invariance, the frequency structure of the PSF at the 
zenith pixel is same as the sidelobe structure of the 
PSF in I evaluated at f Q . The frequency structure of 
the PSF in a bandwidth Sf is identical to the sidelobe 
structure of the PSF (evaluated at f Q ) in the interval 
[l -6fl /(2f ) , l +Sfl /(2f )]. This truncated ver- 
sion of the PSF can be extracted by multiplying the 
PSF with a window centered at l D and having width 
Sfl /f . A Fourier transform of the product of these 
two functions is the convolution of their individual 
Fourier transforms which are the baseline weighting 
function W(x) and the convolution kernel from Equa- 
tion [TT] It may be seen from Equation \TU\ that the 
sinc-function envelope of the convolution kernel has 
a main lobe width equal to 2c/(/ Q A/), which in gen- 
eral is far narrower than the visibility range 2x max . 
The rapid variation in the convolution kernel is a com- 
plex sinusoid in x with frequency l f /c. At any 77, 
the convolution in Equation [TT] is sensitive to a small 
range in W{x) that is centered at x = cq/l with a 
width c/(l a Af) which is the main-lobe width of the 
convolving kernel. For a visibility distribution with 
uniform weight and complete coverage like the one dis- 
cussed, owing to the rapid fluctuations of the complex 
sinusoid, there exists significant response only around 
Xmax- This explains the delta- function like 'PSF con- 
tamination' in rj centered at 77 = x max l a /c (see Fig- 
ure [TJ . Sources from different directions in the sky 
have different values of l Q and hence contribute to fore- 
ground contamination at different fairly localized val- 
ues of rj. The relative magnitude of these impulses is 
determined by the relative sidelobe levels at the re- 
spective directions. Figure [3] shows the delta-function 
like 'PSF contamination' in 77 for point sources entering 
the imaged pixel at zenith through different parts of 
the PSF (different l a ). A rectangular window has been 
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assumed in frequency domain. The finiteness of W(x) 
gives a 1/2 roll-off in sidelobe levels which corresponds 
to a I/77 roll-off in the contamination at different values 
of 77. 
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Fig. 3. — 'PSF contamination' in 77 space due to 1 Jy 
point sources at l a = 0.05,0.1,0.2,0.4, and 0.5. The 
left panel shows the contamination for W{x) with uni- 
form weight and a rectangular frequency bandpass, 
and the right panel shows the contamination for a 
Blackman-Nuttall window function. For a visibility 
distribution with complete coverage, the PSF sidelobes 
are a result of the 'top hat' shape of the weighting func- 
tion W(x), and the convolution in Equation [TT1 has a 
delta-function like impulse located at r\ = x max l /c. 
The magnitudes of such impulses are equal to the 
sidelobe level at the value of l where they are eval- 
uated. The solid black line is an approximate I/77 fit 
to the peaks of the delta-function like impulses. For 
V > Xmaxl/c, the sidelobes of the convolving kernel 
give a response that further falls off as I/77. The 
Blackman-Nuttall window substantially lowers con- 
tamination in the 'EoR window' (77 > x max l /c). 



3.2. Sidelobes of the convolving kernel 

We have used the phrase 'delta-function like' to de- 
scribe the 'PSF contamination' for a source entering 
the zenith pixel from the direction l a . The exact pro- 
file is described by the sine envelope of the convolving 
kernel. Owing to this, the 77 space may be divided 
into two regions. For l a G [01], the delta-function like 
peaks are confined to 77 £ [0 x max /c\. This defines 
a region of relatively high foreground contamination, 
say Hi. However, for 77 > x rnax /c, there still ex- 
ists 'PSF contamination' due to the sidelobes of the 
convolving kernel. The contamination in this regime, 
say 7^2, falls off as I/77 since the sidelobes of a sinc- 
function envelope of the convolving kernel have a 1/x 
roll off. While EoR power spectrum estimation for 
77 G [0 x max /c] will suffer from stronger contamina- 
tion, the level of contamination for 77 > x max /c may 
be reduced by a judicious choice of a bandpass window 



B(f) that gives low sidelobe levels in 7^2, and hence, 
IZ2 defines an 'EoR window' where we may focus our 
efforts to detect EoR. 

In essence, a Fourier synthesis array gives a PSF 
that varies with frequency. If the baseline extent of 
the array is x max , the spectrum of frequency varia- 
tions will have significant energy only in IZi defined 
by 77 G [0 , x max /c]. The instrument will sample 
the ?7-space up to r/ max — 1/23. In 72.2 defined by 
V G [xmax/c , 1/S] the foreground contamination is 
due to the l/77-form sidelobes of the contaminants in 
IZi. These I/77 sidelobes are the result of a 'top hat' 
frequency window. A smooth window function in fre- 
quency like, fo r example, a Blackman-Nuttall window 
(|Nuttall 1119811 ) will substantially reduce the sidelobe 
contamination in 72-2 (by more than three orders of 
magnitude) at the cost of reduced resolution in 77 space 
(see Figure E]) . 

3.3. Instrumental k space 




k ± , Mpc" 1 



Fig. 4. — Figure showing instrumental response in k 
space at redshift of z = 8. The MWA instrument 
parameters have been used for illustration. A band- 
width limit of 30 MHz defines a minimum feii m j n of 
about 0.012 Mpc -1 , a frequency resolution of 40 kHz 
defines a k\\ max of about 9.1 Mpc -1 , an array extent 
of 1 km defines a k± max of about 0.3634 Mpc -1 , and a 
minimum baseline of 7.5 meters gives a k± m i n of about 
3.5 x 10 -4 . The region of high foreground contamina- 
tion IZi extends up to rj — x max /c — 3.33 /isec that 
corresponds to k\\ m of about 1.2 Mpc -1 . The region 



between kii m and k 



II' 



is the 'EoR window' IZo and 



is a region of low foreground contamination. The un- 
filled black circles define the point (k± ma x, k\\ m ) for 



redshifts of 3,4,5,6,7,8,9, and 10 
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The instrumental k space is an effective space to 
represent regions and extent of foreground contamina- 
tion. The instrumental k space is 2 dimensional with 
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the comoving line of sight wavenumber fen and comov- 
ing transverse wavenumber k± as axes. If k x , k y , and 
k z are comoving wavenumbers along three spatial axes 
with origin at some point ab out which the EoR fluctua- 
tions will be estimated, then lMorales fc Hewitt! (|2004) 
have shown that 



k± — \/ k^ + ky, k x 
2ttv 



2ttu 

D M { z y 

2TrH f HI E(z) 



D M {zY 



c(l + zf 



(12) 



where the approximation holds for a small bandwidth 
(small redshift range). Here Jhi = 1420 MHz is the 
rest frequency of 21 cm HI line emission, Hq — 70 km 
sec -1 Mpc -1 is the present value of the Hubble con- 
stant, E{z) = (Om(1 + z a f + n k {l + z Q f + 0\) 5 , 
and Dm{z) is the transverse comoving distance at red- 
shift z. Figure Q] shows a plot of the measurement 
space (in k± — k\ | domain) for a Fourier synthesis array 
with complete visibility coverage. A Fourier synthesis 
array is sensitive to measurements in a region in the k 
space bounded by some maximum and minimum val- 
ues of k± and k\ | . The highest k± mode sampled by 
the instrument is related to the uv extent of the array 
and is given by 



27TX max f J 
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c(l + z)D M ' 



(13) 



while the smallest k± probed by the instrument is re- 
lated to the antenna extent d, and is given by 
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(14) 



The highest value of k\\ sampled by the instrument is 
related to the frequency resolution and is given by 



*ll 



max ~ yjmax 



2nH f H iE{z) 
c(l + z) 2 



3 7]max 



l 

A/' 



(15) 



The lowest value of k\ | sampled by the instrument in a 
snapshot observation is related to the frequency band- 
width and is given by 



min ~ y]min 
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(16) 



The final comment in this section is related to the 
implications of If invariance of the PSF to our un- 
derstanding of the location and level of contamination 
in the instrumental k space. Consider the foreground 
contamination in the neighborhood (61 << 1) of some 
image pixel V due to a point source at some l removed 
from V. Due to the // invariance, the transverse struc- 
ture of contamination around V is identical to the line 
of sight contamination at V. Since ^ = y-, the en- 
ergy in transverse fluctuations around V at wavenum- 
ber &_lo will be the same as the energy in line of sight 



fluctuations at V at wavenumber r\ Q = -j- k± . Using 
Equation ll21 the (fcj_, fcii) pairs with identical contam- 
ination in them are defined by fe|| = k± 2 ^°^^ l a . 
For a uniformly weighted visibility distribution with 
complete coverage, there is significant energy only in 
the transverse mode k± = x max fo/c- This gives us a 
corresponding rj a = l a / f a k± — ^max^o/ C- & result we 
established in Section l3~Tl 



4. Effect of gridding and imaging algorithms 

While Section [3] described the nature of foreground 
contamination along the line of sight and transverse 
wavenumber axes, its results are a good description 
of the contamination in the special case of an ideal 
array, and hence lay the foundations needed to ap- 
preciate issues related to EoR detection with practi- 
cal array geometries and data processing algorithms. 
Real arrays seldom have complete visibility cover- 
age, and the Fourier relationships between visibility 
and image space, and between frequency and line of 
sight wavenumber space are, in practice, evaluated as 
discrete transforms. Computationally efficient trans- 
forms like the Fast Fourier Transform (FFT) are often 
used in Fourier synthesis arrays and almost always 
involve re-gridding the visibility data on a regular uv 
grid. This section discusses the effects of such gridding 
and Fourier transform algorithms used for imaging on 
the nature of the 3 dimensional PSF and the result- 
ing foreground contamination in k space. The first 
generation EoR experiments propose to differentiate 
confusing foregrounds from the EoR signal by attempt- 
ing to synthesize a confusion-free 'EoR window' in the 
line of sight wavenumber dimension. We now describe 
the effect of gridding and imaging algorithms on the 
foreground contamination in rj. 

We have shown in Section[3]that contamination in 77 
space is substantially different in the two regions (see 
Figures [3] and [4| IZi and IZ2 , where IZi extends from 
rj = to 77 = x max /c, and IZ2 (or 'EoR window') ex- 
tends from rj — x max /c to 77 = 1/A/. We have also 
shown that contamination in 7?-2 is due to the side- 
lobes of the convolution kernel B (see Equation ITT]) 
and the delta-function like peaks from sources at var- 
ious directions are confined to 72>i . Consequently, if 
region H2 exists (1/A/ > x max /c), as it should in 
any well-designed EoR instrument that aims to syn- 
thesize a confusion-free 'EoR window', then PSF vari- 
ation in frequency at any direction cosine 1 is oversam- 
pled by the instrument spectral resolution A/. This 
implies that, in the ideal case, there will not be sig- 
nificant channel to channel variation in flux versus 
frequency at any sky pixel. However, processes like 
gridding may induce channel to channel 'jitter' in the 
PSF and potentially contaminate H.2- This stochastic 
component of foreground contamination arising from 
the gridding process may be termed 'gridding contam- 



G 




ination'. We now describe such stochastic effects of 
gridding in terms of channel to channel migration of 
baseline vectors from one uv pixel to another. 

4.1. Baseline migration 



Ch#n 



Ch #n+l 



O Native grid point (irregular grid) 
B Post gridding convolution grid point (regular grid) 
f\ Convolution kernel 

Fig. 5. — Depiction of inter-channel baseline migra- 
tion in one dimension (it). Red squares are grid points. 
Blue circles are baseline vector locations in wavelength 
units. The continuous line represents the post grid- 
ding convolution visibility distribution for a triangular 
convolution kernel. For simplicity, we have assumed 
that both depicted visibility measurements have equal 
weights. The upper panel shows the visibility axis for 
some frequency channel n, and the lower panel shows 
the visibility axis for the adjacent frequency channel 
n + 1. The depicted visibility segment is shown to 
have a fairly sparse u coverage for clarity. Notice how 
the baselines migrate on the visibility grid. While grid 
pixel number 1 and 8 go from being 'filled' to being 
a 'hole' between channel n and n + 1, visibility pixels 
5 and 12 go from being a 'hole' to being 'filled'. Pix- 
els like 2 and 10 continue to be filled but experience a 
'step' change in their weighted visibility value. 

Fourier synthesis arrays have traditionally formed 
image cubes in Imf space by (i) computing the base- 
lines (in wavelength units: u = x.//c) generated by 
a given antenna distribution at every frequency chan- 
nel, (ii) weighting the visibilities to get desirable PSF 
characteristics, (iii) performing gridding convolution 
on the visibilities to re-grid them on a regular uv grid 
at every frequency channel, and (iv) using an efficient 
FFT algorithm to compute a 'dirty' sky map at ev- 
ery frequency channel. Among these signal processing 
algorithms, gridding convolution has not been con- 
sidered in the analytical treatment of Section [2] We 
now describe the influence of gridding convolution on 
foreground contamination in 77 space. 

A baseline vector x (in meter units) will suffer a 
displacement of <5u = x<5//c (in wavelength units) be- 
tween adjacent frequency channels in the uv plane. 
Here Sf is the inter-channel frequency spacing. We 



call this phenomenon 'baseline migration'. This dis- 
placement results in 'step' changes in gridded visibility 
values of pixels within the support of the convolving 
kernel. Many such 'step' changes in visibilities due to 
the displacement of many baseline vectors results in 
a stochastic 'jitter' in the frequency structure of the 
PSF. PSF 'jitter' due to baseline migration may in- 
duce channel to channel fluctuations in flux versus fre- 
quency on the image plane, thereby significantly con- 
taminating the 'EoR window' TZ 2 . This PSF 'jitter' 
due to 'baseline migration' is expected to have high 
contribution (i) from longer baseline vectors due to 
the relatively higher quantum of baseline displacement 
<5u, and (ii) from regions of low uv density. In regions 
of higher uv density a given uv pixel, during grid- 
ding convolution, receives contributions from several 
neighboring baselines. The visibility contribution to a 
given uv pixel then depends on the relative displace- 
ment of the neighboring baselines, rather than their 
actual displacement. On the other hand, consider the 
extreme case of very low uv density as depicted in 
Figure [5] Here the actual inter-channel displacement 
of the baseline vectors induces relatively high 'steps' 
in the post-gridding convolution visibilities. Natural 
weighting results in significantly lower weight assigned 
to regions with (i) longer baseline vectors and (ii) low 
uv density and is expected to result in low PSF 'jitter' 
in frequency. Uniform weighting, on the other hand, 
increases the weighting for longer baselines and is ex- 
pected to result in substantially high PSF 'jitter'. 

4.2. The central uv void 

The angular power spectrum of Extragalactic 
source s has been observe d to have a power-law like 
form ([Blake fc Walll 120021) with substantially greater 
power on small angular scales; the universe is homo- 
geneous and isotropic on large scales. Consequently, 
short baselines are not expected to have substantially 
greater response to Extragalactic continuum sources. 
However, the mean sky brightness will have a response 
at zero spacing (u = 0,v = 0). Energy in various spa- 
tial Fourier modes is seldom perfectly isolated and the 
zero spacing component will leak into short baselines. 
Consequently, channel to channel changes in the visi- 
bility distribution close to u = 0, v — may result in 
'step' changes in flux versus frequency at image pix- 
els thereby contaminating higher values of 77 in 1Z2- 
While Extragalactic sources are expected to have a 
smooth distribution on large angular scales, flux from 
the Galactic synchrotron emission is expected to follow 
a power law with significant spatial distribution power 
at large angular scales. Consequently short baselines 
close to u = 0,v — will have significant response to 
our Galaxy and 'baseline migrations' in this part of 
the uv plane will create significant PSF 'jitter' that 
may contaminate 72.2- 

Fourier synthesis arrays usually have a central l uv 
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void' as the antenna elements have a minimum spac- 
ing due to constraints imposed by their own physical 
size. Antenna cross-talk considerations often result 
in array configurations with a central 'uv void' that 
is several times the antenna size. In such cases, the 
central uv void spans several uv pixels. We mentioned 
in Section H3] that 'baseline migrations' in regions of 
high uv density do not contribute significantly to inter- 
channel PSF 'jitter'. Though short spacing baselines 
lie in a densely populated part of the uv plane, with 
increasing frequency, most baseline vectors are dis- 
placed radially away from any uv grid point close to 
u = 0, v = 0, as opposed to other uv grid points where 
baseline vectors also get displaced towards them. This 
asymmetric neighborhood near short spacing visibil- 
ities results in relatively higher inter-channel 'step' 
changes in the visibilities on short baseline vectors. 
To alleviate this problem, the central uv void may 
be chosen to be larger than that imposed by antenna 
cross-talk constraints. 



4.3. The Chirp Z transform 

We have seen how channel to channel changes in 
the visibility weighting W(u, /) contaminates the line 
of sight wavenumber dimension. Here W(u, /) in- 
cludes the effect of all visibility weighting and tapers 
used, and W(u, /) determines the PSF via the imaging 
transform. If the visibility distribution has complete 
coverage, maintaining invariance of W(u, /) with fre- 
quency is possible through application of a frequency 
independent visibility taper function. In such a case, 
W(u, f) is same at all frequencies and there is no 
frequency dependence in the PSF In practice, there 
are m any holes in the uv grid and lLabropoulos et al. 
(2009) have suggested that a visibility mask may be 
used that rejects the uv pixels that are not filled at 
all (or most) frequencies prior to image formation. 
At one end, we may exploit the complete instru- 
ment data and sensitivity at the cost of increased 
contamination in 7^-2, and at the other end we may 
have a frequency independent PSF with some loss 
in information and sensitivity. However, there ex- 
ists a 'middle path' that assures that no data is re- 
jected while restricting the foreground contamination 
in IZ2 — i mage formation using the Chirp Z Transform, 
or CZT dRabiner. Schafer fc Raderl Il969h . We first 
present a brief description of CZT. 



The Z transform of any finite discrete sequence 
x n , n — 0, 1, 2....N — 1 is defined as 



n=N-l 



x n z 



(17) 



n=0 




2 [2n/(Na)] + arg(A) 

1 [2n/(Na)] + arg(A) 
+ arg(A) 



Z transform is a function of the complex variable 
z. When evaluated at discrete points along the 
unit circle defined by Zk = exp (i2nk/N), where 



(N-l) [2n/(NcO] + arg(A) 
(N-2) [2n/(Na)] + arg(A) 

Fig. 6. — Sketches showing points on the z plane on 
which the generalized Chirp Z transform is computed 
for two cases: Case 1 (left panel) is for A — 1 and W = 
exp (— i2n/N) and Case 2 (right panel) is for \A\ = 1 
and W — exp (—i2ir/ (Na)). In Case 1, the generalized 
CZT reduces to DFT, and for Case 2 the generalized 
CZT reduces to a special case of CZT that will be used 
in the proposed imaging algorithm of Section 14.41 



k = 0, 1, 2....N — 1, Z Transform of the finite discrete 
sequence x n reduces to its Discrete Fourier Transform 
(DFT). The Chirp Z Transform or CZT is another 
special case of the Z transform that is evaluated at 
discrete points on the locus given by zu = AW~ k 
where A and W are complex constants. Such a locus 
describes a spiral in the complex z plane. Notice how 
for A = 1 and W = exp (-i2ir/N) the CZT again 
reduces to a DFT. 

Consider a case where \A\ = 1 and W = exp (—i2ir/ (Na)) 
where a is a real constant. The locus of points on the 
complex z plane where the transform is evaluated now 
lie on a unit circle like in the case of a DFT. However, 
angular spacing between the points on the locus is 
now 2ir/(na) as opposed to 27r/A^ in case of a DFT. 
The first point on the locus is also shifted by the ar- 
gument of the complex constant A. The loci for the 
two cases are shown in Figure [51 Though \A\ = 1, 
W — exp (—i2n/(Na)) is a special case of the gener- 
alized Chirp Z Transform, we refer to it hereinafter as 
Chirp Z Transform or CZT with a chirp factor a. The 
only difference between DFT and CZT is the presence 
of the chirp factor a and an offset due to the pres- 
ence of A. a simply 'stretches' the locus on which 
the transform is evaluated and the CZT of a sequence 
is essentially a DFT of a 'stretched' version of the 
sequence — the stretching factor being a. We now de- 
scribe how visibilities can be gridded in meter units 
and how a chirp factor can be used to compensate for 
the 'stretching' of baseline vectors with frequency. 



4.4. A CZT based imaging algorithm 

As described in Section 14.11 the principal cause 
of channel to channel PSF variations in conventional 
gridding and imaging processes is baseline migration 
due to independent visibility gridding at every channel. 
However, the baseline vector generated by a given an- 
tenna pair moves across the visibility grid with change 
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in frequency in a very predictable manner. The base- 
lines in wavelength units at any frequency / are merely 
the baselines in meter units 'stretched' by a factor a, 
where a = f/c. In principle, we may grid the visibil- 
ities only once in meter units and use the CZT with 
chirp factor a to form the image at every frequency 
channel. Given an image extent, the grid size in meter 
units is chosen to avoid aliasing even at the highest 
frequency in the visibility cube. Gridding the visibili- 
ties only once in meter units eliminates the problem of 
non-linear baseline migrations and consequently ame- 
liorates the problem of contamination in IZ2 . In sum- 
mary, we grid the baselines in meter units, and thereby 
do not let the baselines scale with frequency, and a 
given antenna pair will fall on exactly the same visi- 
bility grid point at all frequencies. We then compute 
the image on an Im grid that is constant over frequency 
and use the chirp factor a in a CZT to compensate for 
an unchanging baseline. Such a CZT based imaging 
transform is given by 



'S'(lmn) 



=N-lj=N-l 

j2n 



exp 



Na(f) 



where i, j are uv (meter units) grid pixel indicies, to, n 
are image grid pixel indicies, function Vi(xjj, /) repre- 
sents the visibilities at frequency / gridded on to the 
common grid, and a(f) is the chirp factor given by 
a(f) = c/f. The transform in Equation 1181 is simi- 
lar to a DFT and may be evaluated by the same ef- 
ficient FFT algorithms with a minor modification to 
the exponential terms. The CZT image formation re- 
sults in reduced contamination in IZ2 and, moreover, 
at a substantially reduced computational cost of grid- 
ding convolution because the gridding convolution at 
all frequency channels is done to a single common grid 
with little additional computational cost for imaging. 

5. The MWA example 

We have so far discussed the nature of foreground 
contamination in r) space for an ideal array with com- 
plete visibility coverage. We have also outlined some of 
the practical considerations applicable to realistic ar- 
rays and in the process proposed a new imaging tech- 
nique using the Chirp Z Transform. In this section we 
demonstrate the relationships and results of Sections 
[21 El and H] using the case of snapshot imaging with 
MWA. MWA consists of 512 antenna elements called 
'tiles'. Each tile consists of a phased array of 16 bow- 
tie antennas. The tiles are distributed to form base- 
lines up to about 1 km. The instantaneous bandwidth 
is about 30 MHz, and we assume a frequency resolu- 
tion of about 40 kHz. MWA will produce snapshot 
images every 8 seconds. The instrumental response 
for this instrument configuration in the k± — fci plane 



is shown in Figure [4] Details of the MWA te lescope 



(2009). 

5.1. Confusing sources and sidelobes 

To estimate ' blending confusion', we u s e the source 
counts given by ISubrahmanvan fc Ekersl f|2002b . The 
threshold flux density Sc above which we expect to 
find one source per beam of full width at half maximum 

" arcmin 

at frequency f GHz is given by 



S c = 40 



3 1.67 



arcmin J GHz 



ft 



-0.67 



juJy. 



(19) 



Sc is a confusion limit for a survey with excellent visi- 
bility coverage and all sources much stronger than Sc 
can, in principle, be identified and subtracted. The 
stochastic response to sources weaker than Sc con- 
tribute to sky pixel rms that is close to Sc- Sc for 
MWA evaluates to 2.3 mJy and is about 3 orders of 
magnitude higher than the expected thermal noise in 
8-second snapshot images: MWA snapshot images will 
be confusion limited. For our simulations we assume 
successful peeling of point sources above 5Sc ~ 10 
mJy. The expected residual rms confusion noise af- 
ter successful subtraction of sources above 10 mJy is 
about 3.6 mJy. Due to the dominance of confusion over 
thermal noise, we will neglect the effects of thermal 
noise hereinafter. However, it must be noted that post 
foreground modeling and successful subtraction, the 
residual image that will be used for EoR power spec- 
trum estimation is expected to be dominated, among 
other residual fluctuations, by thermal noise. 

The threshold to which continuum sources may be 
identified and subtracted also depends on our ability 
to distinguish spillover flux due to sidelobes on the 
sky. We now roughly estimate the expected rms fluc- 
tuation at any image pixel due to the sidelobes on all 
the confusing sources in the sky — a quantity we call 
'sidelobe confusion'. If we assume (i) the sidelobes 
of the PSF to be a zero mean random variable, and 
(ii) the confusion flux in each pixel to be a random 
variable that is uncorrelated with the sidelobe level, 
then the approximate rms sidelobe confusion can be 
estimated fairly easily: Cs = Y^i where Pi is the 
sidelobe level at pixel i, and Si is the primary beam 
weighted true sky flux due to confusing sources at pixel 
i. Since (Pi) = 0, the variance in sidelobe confusion is 
{CD = E t (P?S?) = £ 4 (/f ) (Sf). The far sidelobes 
of MWA are expected to have an rms variation of 
about 0.05% giving (Pf) = 0.25xl0~ 6 , and the RMS 
confusion is about 3.6 mJy giving (Sf) = 13 mJy 2 . 
The sidelobe level assumed here is true for sidelobes 
far from the phase center at which 'sidelobe confusion' 
is being computed. In reality sidelobe levels vary from 
sky pixel to pixel and a constant value has been as- 
sumed here for simplicity. Far sidelobe level has been 
chosen since far sidelobes result in foreground con- 
tamination at higher values of r\ where we expect to 
detect the EoR with minimal contamination, and this 
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calculation will present the worst case. While Pi is as- 
sumed to be independent of i, Si is assumed to be the 
confusion noise weighted by the primary beam gain, 
and the summation is performed using the expected 
primary beam of the MWA primary antenna element 
(Full width to first null of 60°) to evaluate 'sidelobe 
confusion'. The resulting 'sidelobe confusion' from 
the entire sky is (Cg) « 1.9 mJy. The similarity in 
the values of 'blending confusion' and 'sidelobe con- 
fusion' is an interesting consequence of the similarity 
between the number of resolution elements in the im- 
age (weighted by the primary beam) and inverse of 
the sidelobe level of the synthesized beam. This im- 
plies that after subtracting bright point sources from 
MWA snapshot images, the resulting residual image 
has statistical fluctuations with roughly equal contri- 
butions from the stochastic response to the Extra- 
galactic sources and the sidelobe levels. 

The EoR brightness fluctuations are expected to 
have an rms of about 25 mK on scales comparable to 
the MWA 5 arcmin resolution. This corresponds to an 
rms flux of about 36 fiJy. Hence, apart from confusing 
foreground sources, 'sidelobe confusion' can confuse 
any attempt to detect the EoR signal through a simple 
measurement of image plane variance. That said, flux 
from confusing sources and PSF sidelobes have struc- 
ture in the frequency domain which can be exploited 
to detect the EoR fluctuations despite a high rms fluc- 
tuation in the image plane. While spectral structure of 
confusion Extragalactic sources follow relatively sim- 
ple power laws, spectral structure of instrumental side- 
lobe confusion is more complex as described in Sections 
[2j O and SJ We now quantitatively estimate the fore- 
ground contamination in frequency space (and eventu- 
ally in rj space) for the case of snapshot imaging with 
MWA. 

5.2. A snapshot observation 

This subsection describes a simulation to estimate 
the foreground contamination in r\ space for the case of 
snapshot imaging of Extragalactic point sources using 
MWA. The simulation assumes a frequency invari- 
ant sky model and assumes that all sources brighter 
than 10 mJy have been successfully subtracted. The 
sky model is then interpolated into an Im grid with 
grid size at least 2 times finer than the MWA instru- 
ment resolution. The grid extent is covers the interval 
[—0.5 0.5] in direction cosines that corresponds to the 
MWA primary beam full width at first null, which is 
60° at 150 MHz. Gridding convolution of visibilities 
is performed using a triangular kernel like the one 
shown in Figure [3] A simple kernel was chosen purely 
to illustrate the nature of contamination in rj space. 
At every frequency channel, three separate PSFs are 
generated using natural weighting, uniform weighting, 
and uniform weighting with a Gaussian taper (here- 
inafter referred to simply as Gaussian weighting) . The 



simulated sky in Im is then multiplied with the PSF 
at zenith and the result is summed. This accumulated 
value at every frequency channel gives the dirty image 
flux density at the zenith pixel. This zenith pixel flux 
versus frequency is used to evaluate the spillover of 
foreground contamination by Fourier transforming to 
r\ domain. The simulation was performed for three 
cases each involving a different gridding and imaging 
algorithm; we discuss these below. 

Case 1 employed independent gridding convolution 
at every frequency channel, and the dirty image was 
formed using an efficient FFT algorithm. Figure [7] 
shows the flux density versus frequency and contami- 
nation power in rj for this case. We have shown earlier 
that for rj < x max /c in 7£i, the exact structure of 
spillover depends on the nature of relative weightings 
given to visibilities. Weighting here includes the ef- 
fect of uv plane distribution of baselines, user defined 
weighting and any tapers applied. For the case of 
uniformly weighted visibilities with complete uv cov- 
erage, the rj spillover from a point source is expected 
to be localized since the convolution in Equation QT] 
has a significant contribution only due to the disconti- 
nuity in the weighting function at |x| = x max . In case 
of MWA snapshot visibilities, the weighting function 
is not 'smooth' for |x| < x max and the convolution 
in Equation Qj] has a significant contribution from 
all |x| < x max . Though this results in a distributed 
spillover, a significant part of the spillover energy is 
still restricted to TZx, and IZ2 will, in the absence of 
baseline migrations, be contaminated only due to the 
sidelobes of the convolving kernel from Equation [11] 
as discussed earlier. This is the 'PSF contamination' 
that is expected to persist in TZ 2 - Any additional chan- 
nel to channel variations in the PSF will contaminate 
IZ2 beyond the 'PSF contamination'. Such channel to 
channel variations are evident in the upper panel of 
Figure [7] in the form of 'jitter', especially for the cases 
of uniform weighting and Gaussian weighting. This 
'jitter' was found to be due to inter-channel baseline 
migrations as discussed in Section^] To quantitatively 
estimate the amount of 'gridding contamination' due 
to baseline migration, the Blackman-Nuttall window 
with very low sidelobes was used as a frequency win- 
dow function B(f). As described in Section I3T21 such 
a windowing function will significantly suppress 'PSF 
contamination' in IZ2 (r/c > 500 meter), and the dom- 
inant residual contamination in IZ2 will be due to the 
inter-channel 'jitter' arising from the gridding process. 

It is interesting to note that baseline migrations 
have induced relatively higher jitter for the cases of 
uniform weighting and Gaussian weighting as opposed 
to natural weighting. As discussed in Section 14.11 this 
is due to two reasons. First, for an inter-channel fre- 
quency spacing of Sf, a baseline vector x (in meter 
units) suffers an inter-channel migration on the uv 
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plane of length xSf/c (in wavelength units). Second, 
and more importantly, a relatively higher amount of 
inter-channel 'jitter' due to baseline migrations arises 
from regions of lower uv density (longer baseline vec- 
tors) compared to regions of higher uv density (smaller 
baseline vectors). Since uniform weighting, and in this 
case Gaussian weighting, place higher emphasis on 
longer baseline vectors as compared to natural weight- 
ing, the amount of 'jitter' in Figure[7]goes in decreasing 
order of magnitude from uniform weighting to Gaus- 
sian weighting to natural weighting. Earth rotation 
synthesis and drift scan-synthesis will increase the uv 
density and is expected to ameliorate the 'jitter'. 

Case 2 involves image formation with the CZT al- 
gorithm with a similar image resolution at every fre- 
quency. As described in Section 14.41 image formation 
with the CZT algorithm is expected to lower 'gridding 
contamination' in IZ2 arising from inter-channel jitter 
due to 'baseline migrations'. Figure [3] shows the line 
of sight contamination for the case of image formation 
using the CZT algorithm. The 'jitter' in flux density 
versus frequency (upper plot) has significantly reduced 
due to the absence of baseline migrations inherent to 
independent gridding at every frequency channel — an 
expected outcome of gridding the visibilities in me- 
ter units. Similar image resolution at every frequency 
channel was achieved by scaling the Gaussian taper to 
have the same form in wavelength units, and restrict- 
ing the longest baseline vector that was used in the 
gridding and imaging process to the same value (in 
wavelength units) at every frequency channel. This 
resulted in inter-channel 'steps' in the weighting func- 
tion at the far end of the uv plane that resulted in 
a small but significant amount of 'jitter' and hence a 
contamination in 1^2- Nevertheless, the contamina- 
tion in H2 for this case is reduced by a factor of 1.7, 
6.4, and 21.7 for natural weighting, uniform weight- 
ing, and Gaussian weighting respectively. Relative 
contamination for natural weighting has not improved 
significantly as regions in the uv plane that contribute 
most to 'jitter' are down-weighted by natural weight- 
ing. 

Case 3 involves image formation with the CZT al- 
gorithm with a frequency dependent image resolution. 
In this case, the set of baselines that are used in the 
gridding and imaging routine are identical at all fre- 
quencies. It may be noted that despite a frequency 
dependent image resolution, the CZT algorithm com- 
putes the sky flux density on the same Im grid at all 
frequencies. The results of this simulation are shown in 
Figure [9] The contamination in 7?. 2 has dropped con- 
siderably owing to the elimination of channel to chan- 
nel 'jitter' in the PSF by this method. The extremely 
low sidelobes of the Blackman-Nuttall window and 
the absence of 'baseline migration' in the CZT based 
imaging algorithm we have proposed herein result in 



extremely low levels of contamination in H.2- Quanti- 
tatively, the contamination in IZ2 for this simulation 
as compared to Case 1 has reduced by more than three 
orders of magnitude for all types of weighting schemes. 

It is important to note that the amount of inter- 
channel 'jitter' is largely dependent on the visibility 
distribution and in particular on the filling fraction 
in different parts of the uv plane. Therefore, the rel- 
ative merits of the three different weighting schemes 
discussed will depend on the array configuration and 
observing strategy. It is however expected that the 
contamination in TZ.2 will be orders of magnitude lower 
if gridding and imaging are done using the proposed 
CZT algorithm along with a good frequency window. 
The simulations and results in this section are not ex- 
haustive, and are included as a means to appreciate 
the various factors that influence the frequency depen- 
dence of the PSF and illustrate the potential reduction 
in foreground contamination. 

6. Conclusions and future work 

Redshifted 21 cm tomography using Fourier syn- 
thesis arrays has emerged as a promising tool for EoR 
studies. Mitigating the confusion effects of Galactic 
and Extragalactic foregrounds is a major challenge to 
such EoR experiments. Through this work we have 
furthered understanding of contamination arising from 
Extragalactic continuum foreground and instrumental 
effects using a common framework. In particular, we 
have derived analytical expressions for contamination 
along frequency, or equivalently line of sight (LOS) 
dimension, which relates the frequency-axis contami- 
nation in the image cube to the visibility weighting in 
Fourier space. In doing so, we have cast the problem 
in Fourier space where the telescope measurements 
indeed lie, thus enabling us to draw conclusions that 
are directly applicable to instrument design and data 
processing. 

We identified two major sources of foreground 
contamination — 'PSF contamination' and 'gridding 
contamination', and described their structure in LOS 
wavenumber space. The 'PSF contamination' due 
to any point source is localized in LOS wavenum- 
ber space around 77 = x ma xl/c where I is the direc- 
tion cosine of the source with respect to the imaged 
pixel, and x max is the maximum baseline length in 
meters. Consequently, when imaging with contin- 
uum confusion that is offset in direction cosine upto 
a maximum of l ma x from an imaged pixel, most of 
the 'PSF contamination' will be confined to a regime 
TZ,\ defined by 77 <E [0 x max l max /c\. It is useful to 
note that there does exist a regime IZ2 defined by 
V 6 [xmaxlmax / c 1/A/], where Af is the frequency 
resolution, where finite 'PSF contamination' persists. 
This contamination in 7I2 may be suppressed by three 
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orders of magnitude by judicious choice of a window 
function in frequency. For this reason, IZ2 is an 'EoR 
window' where we may focus our efforts to detect EoR. 

We have shown that gridding and imaging routines 
result in additional 'gridding contamination' in LOS 
wavenumber space. In particular, we have shown that 
independent gridding at every frequency channel may 
lead to a stochastic inter-channel 'jitter' in the PSF 
that potentially contaminates the 'EoR window' 7^2- 
To ameliorate problems associated with inter-channel 
'jitter' due to gridding, we proposed an alternative 
gridding and imaging algorithm based on the Chirp 
Z Transform (CZT). In this algorithm the visibilities 
are gridded in meter units rather than wavelength 
units, thereby eliminating 'gridding contamination'. 
The CZT compensates for the stretching of baselines 
(in wavelength units) with frequency by introducing a 
'chirp' term — a scaling factor by which the frequency 
of the Fourier sinusoids are multiplied. We finally 
demonstrated the localization 'PSF contamination' 
and the elimination of 'gridding contamination' in the 
CZT based imaging algorithm using simulations of 
imaging with MWA. 

Our ongoing and future work is directed towards us- 
ing the concepts and ideas in this paper in an all sky 
simulation of MWA images that incorporates different 
gridding convolution kernels, Earth rotation and drift- 
scan synthesis, effects of array non-coplanarity, and 
the Galactic synchrotron emission. We are currently 
developing specific simulation pipelines that will ac- 
complish the all sky simulations with the CZT based 
imaging algorithm. Such simulations will quantita- 
tively estimate not only the merit of the CZT based 
imaging algorithm for EoR detection, but also the ex- 
pected contamination from foreground continuum in 
future MWA data products. 
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Fig. 7. — Foreground contamination along the line of sight dimension in the zenith pixel for the case of snapshot 
imaging with MWA. In generating these representations of foreground contamination, we populate sources in the 
sky out to |l| mox = 0.5. The zenith pixel flux is plotted against frequency (upper panel) and rj (lower panel). Notice 
the jitter in the upper panel for the case of uniform weighting and Gaussian weighting as opposed to the case of 
natural weighting. This jitter results in significant 'gridding contamination' at higher values of r\ as seen in the lower 
panel. Also shown in the lower panel is radial uv density as a function of |x||l| ma2; for the case of natural, uniform 



and Gaussian weighting. The boundary between 1Z\ and IZ2 lies at 77c = |x| max |l|. 
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Fig. 8. — Foreground contamination along the line of sight dimension in the zenith pixel for the case of snapshot 
imaging with MWA using the Chirp Z Transform (CZT). Notice how the jitter evident in Figure [7] has significantly 
reduced in the zenith pixel flux density versus frequency plot (upper panel). Consequently, foreground contamination 
at higher values of r] (lower plot) is significantly lower as compared to that in Figure[7] A finite amount of undesirable 
'jitter' still persists since a different set of baselines vectors were used in the gridding and imaging process at every 
frequency channel so as to keep the final image resolution constant at all frequencies. 
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Fig. 9. — Foreground contamination along the line of sight dimension in the zenith pixel for the case of snapshot 
imaging using the Chirp Z Transform (CZT) with frequency dependent image resolution. Notice how there is 
negligible 'jitter' in the flux density versus frequency (upper panel). The 'gridding contamination' at higher values 
of f] (lower panel) in IZ2 is now wholly owing to 'PSF contamination'. 
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